{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# [HW1] Problem 2: A Simple Classification Approach\n",
    "\n",
    "Import necessary Python packages.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "np.random.seed(0)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Part (d)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def bivariate_normal(X, Y, sigmax=1.0, sigmay=1.0,\n",
    "                     mux=0.0, muy=0.0, sigmaxy=0.0):\n",
    "    \"\"\"\n",
    "    Bivariate Gaussian distribution for equal shape *X*, *Y*.\n",
    "    See `bivariate normal\n",
    "    <http://mathworld.wolfram.com/BivariateNormalDistribution.html>`_\n",
    "    at mathworld.\n",
    "    \"\"\"\n",
    "    Xmu = X-mux\n",
    "    Ymu = Y-muy\n",
    "\n",
    "    rho = sigmaxy/(sigmax*sigmay)\n",
    "    z = Xmu**2/sigmax**2 + Ymu**2/sigmay**2 - 2*rho*Xmu*Ymu/(sigmax*sigmay)\n",
    "    denom = 2*np.pi*sigmax*sigmay*np.sqrt(1-rho**2)\n",
    "    return np.exp(-z/(2*(1-rho**2))) / denom\n",
    "\n",
    "def generate_data(n):\n",
    "    \"\"\"\n",
    "    This function generates data of size n.\n",
    "    \"\"\"\n",
    "    #TODO implement this\n",
    "    ### start 1 ###\n",
    "\n",
    "    ### end 1 ###\n",
    "    return (X,y)\n",
    "\n",
    "def tikhonov_regression(X,Y,Sigma):\n",
    "    \"\"\"\n",
    "    This function computes w based on the formula of tikhonov_regression.\n",
    "    \"\"\"\n",
    "    #TODO implement this\n",
    "    ### start 2 ###\n",
    "\n",
    "    ### end 2 ###\n",
    "    return w\n",
    "\n",
    "def compute_mean_var(X,y,Sigma):\n",
    "    \"\"\"\n",
    "    This function computes the mean and variance of the posterior\n",
    "    \"\"\"\n",
    "    #TODO implement this\n",
    "    ### start 3 ###\n",
    "\n",
    "    ### end 3 ###\n",
    "    return mux,muy,sigmax,sigmay,sigmaxy\n",
    "\n",
    "Sigmas = [np.array([[1,0],[0,1]]), np.array([[1,0.25],[0.25,1]]),\n",
    "          np.array([[1,0.9],[0.9,1]]), np.array([[1,-0.25],[-0.25,1]]),\n",
    "          np.array([[1,-0.9],[-0.9,1]]), np.array([[0.1,0],[0,0.1]])]\n",
    "names = [str(i) for i in range(1,6+1)]\n",
    "\n",
    "for num_data in [5,50,500]:\n",
    "    X,Y = generate_data(num_data)\n",
    "    for i,Sigma in enumerate(Sigmas):\n",
    "\n",
    "        # TODO compute the mean and covariance of posterior\n",
    "        # in the style of ``mux,muy,sigmax,sigmay,sigmaxy = ''\n",
    "        ### start 4 ###\n",
    "\n",
    "        ### end 4 ###\n",
    "\n",
    "        x = np.arange(0.5, 1.5, 0.01)\n",
    "        y = np.arange(0.5, 1.5, 0.01)\n",
    "        X_grid, Y_grid = np.meshgrid(x, y)\n",
    "\n",
    "        # TODO Generate the function values of bivariate normal\n",
    "        # in the style of ``Z = ''\n",
    "        ### start 5 ###\n",
    "\n",
    "        ### end 5 ###\n",
    "\n",
    "        # plot\n",
    "        plt.figure(figsize=(10,10))\n",
    "        CS = plt.contour(X_grid, Y_grid, Z,\n",
    "                         levels = np.concatenate([np.arange(0,0.05,0.01),np.arange(0.05,1,0.05)]))\n",
    "        plt.clabel(CS, inline=1, fontsize=10)\n",
    "        plt.xlabel('X')\n",
    "        plt.ylabel('Y')\n",
    "        plt.title('Sigma'+ names[i] + ' with num_data = {}'.format(num_data))\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Part (e)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "w = [1.0,1.0]\n",
    "n_test = 100\n",
    "n_trains = np.arange(5,205,5)\n",
    "n_trails = 500\n",
    "\n",
    "Sigmas = [np.array([[1,0],[0,1]]), np.array([[1,0.25],[0.25,1]]),\n",
    "          np.array([[1,0.9],[0.9,1]]), np.array([[1,-0.25],[-0.25,1]]),\n",
    "          np.array([[1,-0.9],[-0.9,1]]), np.array([[0.1,0],[0,0.1]])]\n",
    "names = ['Sigma{}'.format(i+1) for i in range(6)]\n",
    "\n",
    "\n",
    "def compute_mse(X,Y, w):\n",
    "    \"\"\"\n",
    "    This function computes MSE given data and estimated w.\n",
    "    \"\"\"\n",
    "    #TODO implement this\n",
    "    ### start e1 ###\n",
    "\n",
    "    ### end e1 ###\n",
    "    return mse\n",
    "\n",
    "def compute_theoretical_mse(w):\n",
    "    \"\"\"\n",
    "    This function computes theoretical MSE given estimated w.\n",
    "    \"\"\"\n",
    "    #TODO implement this\n",
    "    ### start e2 ###\n",
    "\n",
    "    ### end e2 ###\n",
    "    return theoretical_mse\n",
    "\n",
    "# Generate Test Data.\n",
    "X_test, y_test = generate_data(n_test)\n",
    "\n",
    "mses = np.zeros((len(Sigmas), len(n_trains), n_trails))\n",
    "\n",
    "theoretical_mses = np.zeros((len(Sigmas), len(n_trains), n_trails))\n",
    "\n",
    "for seed in range(n_trails):\n",
    "    np.random.seed(seed)\n",
    "    for i,Sigma in enumerate(Sigmas):\n",
    "        for j,n_train in enumerate(n_trains):\n",
    "            #TODO implement the mses and theoretical_mses\n",
    "            ### start e3 ###\n",
    "\n",
    "            ### end e3 ###\n",
    "\n",
    "# Plot\n",
    "plt.figure()\n",
    "for i,_ in enumerate(Sigmas):\n",
    "    plt.plot(n_trains, np.mean(mses[i],axis = -1),label = names[i])\n",
    "plt.xlabel('Number of data')\n",
    "plt.ylabel('MSE on Test Data')\n",
    "plt.legend()\n",
    "plt.savefig('MSE.png')\n",
    "\n",
    "plt.figure()\n",
    "for i,_ in enumerate(Sigmas):\n",
    "    plt.plot(n_trains, np.mean(theoretical_mses[i],axis = -1),label = names[i])\n",
    "plt.xlabel('Number of data')\n",
    "plt.ylabel('MSE on Test Data')\n",
    "plt.legend()\n",
    "plt.savefig('theoretical_MSE.png')\n",
    "\n",
    "\n",
    "plt.figure()\n",
    "for i,_ in enumerate(Sigmas):\n",
    "    plt.loglog(n_trains, np.mean(theoretical_mses[i]-1,axis = -1),label = names[i])\n",
    "plt.xlabel('Number of data')\n",
    "plt.ylabel('MSE on Test Data')\n",
    "plt.legend()\n",
    "plt.savefig('log_theoretical_MSE.png')\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
